%% Figure of DMA's in US Map
clear
clc

%Continental US only?
continental = 1;
 
%Read shape file
S = shaperead('dma.shp')
mapshow(S, 'DisplayType', 'polygon','FaceColor', [0.75 0.75 0.75])

%Choose axes
xlim([-180,-60])
ylim([15,75])
set(gca,'YTickLabel',[],'XTickLabel',[],'visible','off');

if continental==1
xlim([-130,-65])
ylim([20,55])    
end

hold on

% set path here
cd('replication');

%DMA's in sample
T = readtable('data/DMAs_in_sample.csv');
T = table2array(T(:,2));
Names_DMAs = cell(size(S,1),1); %we will store names of DMA's in our sample too.

%Locate where DMA's in sample are located in shape file.
for i = 1:size(S,1)

aux = intersect(S(i).DMA,T); %Is that DMA in our sample?

if isempty(aux)==0    
%Patch our DMA's to a different colour
mapshow(S(i),'DisplayType','polygon','FaceColor','red');
Names_DMAs{i} = S(i).NAME;
end

end

aux=find(~cellfun(@isempty,Names_DMAs))
Names_DMAs = sort(Names_DMAs(aux));
